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SUMMARY 


To advance the computational capability for combustor research, a three- 
dimensional partially elliptic flow computer program is developed from a con- 
cept originated by Spalding. Without requiring three-dimensional computer 
storage locations for all flow variables, the partiadly elliptic flow program 
is capable of predicting three-dimensional combustor flow fields with large 
downstream effects. This program requires only slight increase of computer 
storage over the parabolic flow program from which it was developed. 

A finite-difference formulation for a three-dimensional fully elliptic 
ti.rbuleni- I’eacting flow field is derived. Because of the negligible diffusion 
effects in the main flow direction in a supersonic combustor, the set of finite- 
difference equations can be reduced to a partially elliptic form. Only the 
pressure field is governed by an elliptic equation and requires a three- 
dimensional storage location ir numerical computations. All other dependent 
variables are governed by parabolic equations similar to those in the parabolic 
flow program and require only two-dimensional storage locations. A numerical 
procedure which combines the previously used marching-integration scheme and 
an iterative scheme for solving the elliptic pressure field is adopted in the 
present program. 

Because of a lack of suitable experimental measiirements , calciilations 
were perfoimed to compare with corresponding calculations made using the 



parabolic flow program. Comparisons show that physically meaningful differences 
are predicted. These differences are especially significant in calculations 
of combustion flow fields. Finally, capabilities and limitations of the 
present program are discussed and future extensions are suggested. 
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SYMBOLS 




D^,D^,D^ 


E 

f 

F 

H 

k 

M 


P 

Pr 

S 

T 


u,v,w 

x,y,z 


r 

6 

e 

K 


expression, eq. (3) 

empirical constants associated with the turbulence model 

coefficients of pressure gradients in the x,y,z directions, 
respectively 

empirical constant, eq. (7) 

mass fraction 
expression, eq. (12) 
total enthalpy 
turbulence kinetic energy 
Mach number 
pressure 
Prandtl nijraber 
source term, eq. (l) 
temperature 

velocity components in the x,y,z direction, respectively 

rectangular coordinates 

exchange coefficeint 

distance between two node points 

turbulence dissipation, energy rate 

empirical constant, eq. (t) 
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u 

p 

T 

(P 

Subscripts : 
eff 


e ,w,s ,n,u,d 


E,W,-S,N,U,D 


0 


P 

t 

TH 


w 


viscosity 
density 
shear stress 

general dependent variable, eq. (l). 
effective 

east, west, south, north, upstream, downstream faces, 
respectively. 

east, west, south, north, upstream, downstream points, 
respectively 

water vapor 

jet 

laminar 

arbitrary node point 
turbulent 
total hydrogen 
wall 


4) general dependent variable 

Superscripts : 

+ nondimensional variable, eq. (T) 

* guessed or approximate quantity 

' correction variable. 
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I . INTRODUCTION 


A computational approach to the prediction of three-dimensional turbulent 
reacting flow fields in supersonic combustors along with experimentation has 
been adopted at NASA Langley Research Center for the scramjet combustor develop- 
ment. Such a computational approach is intended to supplement the traditional, 
expensive, cut-and-try experimental approach and to help overcome some of the 
experimental difficulties associated with the scramjet combustor development. 

In the computational approach, the computer storage and computing time 
are always prime concerns in developing a computer program. However, by 
recognizing special characteristics of the flow field, mathematical and physical 
simplifications can often be made, so that useful engineering results can be 
obtained economically within the existing computer capabilities. One such 
computer program is the SHIP^ (Supersonic Hydrogen Injection Program) computer 
program, which has recently been evaluated, improved, and employed for super- 
sonic combustor design at NASA Langley Research Center. ’ The SHIP program 
is based on a set of parabolic three-dimensional flow equations simplified from 
the full Navier-Stokes equations together with equations for turbulence kinetic 
energy, dissipation rate, and species concentrations. Because of the parabolic 
flow simplification, all flow variables require only two-dimensional computer 
storage and a marching integration procedure can be used in the main flow 
direction. The main limitations for combustor development are inability to 
predict recirculation flow fields which occur adjacent to fuel injectors and 
in separated regions near walls, and to predict the downstream feedback effects 
due to the presence of embedded subsonic flow fields. 
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The aforementioned limitations of the SHIP program coiold he overcome hy 
developing a three-dimensional, fully elliptic flow computer program for the 
entire comhustor. Such a three-dimensional computer program would require a 
three-dimensional computer storage location for each of more than twenty flow 
variables. In addition to other mathematical and physical complications, 
such a fully elliptic approach would be impractical and uneconomical using 
present computers. Since recirculations are usually localized, they may be 
either modeled or computed separately from the main combustor flow field. Flow 
separation close to a wall is associated with energy and total! pressure losses 
in the combustor; such losses should be minimized in a practical design. Since 
a predominant flow direction is present in the main flow field of a supersonic 
combustor, (along which the diffusions of mass, momenta, energy, species, 
etc. can be neglected), the usual parabolic approximation can be applied to 
simplify the flow equations. However, embedded subsonic combustion regions can 
be present within the supersonic stream; thus, downstream effects may be felt 
upstream by means of pressure disturbance propagated through the subsonic fields. 
The present partially elliptic flow computer program is developed to take such 
effects into consideration. 

The concept of the numerical partially elliptic approach was originated 
by Spalding. In view of the approximately parabolic nature of the flow field, 
all flow variables except the pressure are stored only in two-dimensional lo- 
cations in the computer. Only the pressure, which is governed by an elliptic 
equation, requires three-dimensional storage. Thus, it is possible to adopt 
an iterative marching-integration procedure whereby several iterations of the 
flow field in the predominant flow direction are made. Each iteration uses an 
improved estimate of the pressure field obtained from the previous iteration. 
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Iterations are continued until a convergent solution is obtained. Such a 
nijmerical procedure, called a partially parabolic procedure by Spalding, was 
used to compute an incompressible turbulent flow in a two-dimensional curved 
duct. The upstream propagation of pressure effects, which cannot be calculated 
by a parabolic procedure, was demonstrated. The computational advantage of 
such a procedure over a fully elliptic procedure is very great and was illus- 
trated in ref. 5* 

The present partially elliptic flow computer program is developed for a 
three-dimensional, mixed subsonic-supersonic, turbulent, reacting flow in a 
rectangular parallelepiped. The governing differential equations and their 
finite-difference formulations are presented in Section II for a three-dimen- 
sional fully elliptic flow field. The necessary boundary conditions for three 
different boundaries are specified. A method of solving these difference 
equations is briefly discussed. In section III, the finite-difference equa- 
tions are reduced to a set of partially elliptic equations. The differences 
between the partially elliptic and the parabolic formulations are discussed, 
and the advantage in computer storage of a partially elliptic program over a 
fully elliptic program is pointed out. In section IV, the numerical computation 
procedures of the present partially-elliptic program are described. In section 
V, numerical applications of the program are performed. Differences between 
the calculations of the present program and of the parabolic program from 
which it was developed are presented. Finally, the capabilities and limitations 
of the present program are discussed and future extensions are suggested. 
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II. FORMULATIONS FOR A FULLY ELLIPTIC FLOW FIELD 


In this section, a finite-difference formulation for a three-dimensional 
fully elliptic flow field is derived. Expressions of boundary conditions for 
three different boundaries are formulated. These formulations could be used 
to develop a three-dimensional fully elliptic flow computer program. For the 
present purpose, however, they are reduced for developing a three-dimensional 
partially elliptic flow computer program. 

To be consistent with the geometry of the present configuration of the 
supersonic combustor, and also with the existing parabolic flow computer pro- 
gram, the formulation is developed based on an Eulerian formulation in a rec- 
tangular coordinate system (x,y,z) with the z-axis coinciding with the main 
flow direction. The flow field is governed approximately by the Navier-Stokes 
equations together with the equation of state and the species equations. To 
account for turbulence effects the laminar exchange coefficients are replaced 
by the corresponding effective exchange coefficients. The turbulence is 
described by the "k-e" two-equation turbulence model. ^ Equilibrium chemical 
reactions are assumed. This assumption is not essential to the problem develop- 
ment", finite-rate reactions with or without unmixedness can be incorporated 
without any difficulties. Thus the governing equations can be written in the 
following general foim: 
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where 4> is a general dependent variable, and is a general exchange co- 
efficient. When 0 = 1, u, v, w, H, f, k or e, equation (l) corresponds, 
respectively, to the continuity, three components of momentum, energy, con- 
centration, turbulence kinetic energy, or turbulence energy dissipation rate 
equation. On the left-hand side of equation (l), the first term represents 
the sum of convection and diffusion in the x-direction, and the second and 

third terms are those in the y and z directions , respectively. The term S , 

<P 

on the right-hand side is called the source term which includes all other 
terms in the differential equation corresponding to (|). The appropriate ex- 
change coefficients F. and source terms for each variable (b are listed in 

<P 4> 

Table I. The general effective exchange coefficient is composed of two parts, 
a turbulent and laminar; i.e. 


r = eff 
(p Pr 


Pr 

eff,(p t,4) 


Pr 


,({) 


where Pr^ Pr^^ , and are, respectively, the turbulent and laminar 

Prandtl numbers and viscosities. The values of s and the constants 
and in Table I are usually determined empirically (for example, ref. l). 

The mmierical formulation is based on a finite-difference form of equa- 
tion (l). A "staggered" grid system (refs. 6, 7 ) is used. An arbitrary node 
point P is surrounded by six neighboring node points denoted by E, W, S, N, U 
and D for east, west, south, north, upstream, and downstream, respectively, 
in three directions as shown in figure 1. The rectangular control volume of P 
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is chosen such that the faces of the control volume are located at the mid- 
points e, w, s, n, u and d between P and its respective neighboring nodes. 

All dependent variables except velocity components are stored at node points , 
velocity components are stored at mid-points between nodes. The projections 
of a typical control volume in the xy and yz planes and the locations for 
velocity components stored are shown in figure 2. 

By integrating equation (l) over the control volumes a set of difference 
equations can be obtained. For example, the net flux of convection and diffusion 
of the general variable (|) in the x-direction through the control volume can be 
written as 


^ - -Op' 


( 2 ) 


with 


= (Ay)(Az) 
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r<t>g (pu)j, 

■ 2 


= (Ay)(Az) 

*vr 
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Ji. 
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(3a) 


(3b) 


The fluxes of convection and diffusion in the y and z directions can be ex- 
pressed in a similar manner. When the source term S, is linearized, the volume 

<P 

integration of the source term becomes 
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S. (Ax)(Ay)(Az) = . + S cj) 

(j) u,(J) p,9 


ih) 


The resultant finite-difference equation is expressed in the following form: 


A d) 
p P 


d)„ + d) + A^ (J) + A^4> 
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d) 
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A^ (p 
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u,(j) 


( 5 ) 


with 


A 

P 


A^ + A^ + A^ + A^ + A^ + A^ 
e w s n d u 



( 6 ) 


The present formulation is made for a flow field in a rectang-ular parallele- 
piped. To solve the set of difference equations of form (5)? the "boundary con- 
ditions for all variables <}) must be specified at all bour.daries. In general, 
each of the six surfaces of the parallelepiped may be either a surface with uni- 
form flow (e.g. freest re am ) , a symmetry boundary, or a solid wall. At the 
boundary with uniform flow or of symmetry, the uniform or symr.etiy boundary 
conditions are specified respectively. For wall boundaries, the variations of 
flow variables are very steep close to walls. To avoid the use of extremely fine 

Q 

grids near walls, the wall function method is used. Since the present turbulence 
model is valid for fully turbulent flows, all node points (except those at walls) 
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must be located in fully turbulent flov. Hence, the near wall points must be 
far enough from the wall to ensure that the local Reynolds number is 
sufficiently large. The flow between the near wall points and the wall is 
assimied to obey the law-of-the-wall. For example, the velocity component w 
is expressed as 


+ 


w 


1 + 

- Jin E y 
K ^ 


(T) 


where w"*"=w/(t/p)^'*^^ and y'*'=y(Tp or pci'^^k^^^6/Uo at the near wall points 

in the present turbulence model. The empirical constants K, E and are given 
values of 0.ii2, 9-0 and 0.09 respectively; 6 is the distance of the near point 
to the wall. Thus the wall boundary condition for each variable, except the 
turbulence dissipation rate e, can be formulated in terms of the flux of the 
variable (at the near wall point) with an appropriate exchange coefficient 
derived from the law-of-the-wall. The diffusion flux of e to the wall is 
difficult to express. Instead of using the flux and an appropriate exchange 
coefficient, the dissipation rate itself is specified at the near wall point. 


e 


near waJ.l 




(8) 


A summary of the wall boundary conditions for all dependent variables 0 at 
the near wall points is presented in Table II. 

Several computation procedures are available for solving these three- 
dimensional finite-difference equations. A relaxation procedure used by 
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This 
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Spalding et al has been found to be efficient and successful. ’ 

procedure involves a guess and correct scheme for the coupling of velocity with 
11 

pressure , a "hybrid” scheme for modifying the coefficients of the difference 
12 

equations , and a triple sweep line-by-line iterative scheme using a standard 
tridiagonal matrix algorithm in the x, y and z directions. Details are given 
in references 7» 9 and 10. It should be pointed oiit that the solution of the 
finite -difference equation (5) at an arbitrary point requires solutions at six 
neighboring points. Therefore, three-dimensional storage locations are required 
for all variables. Such a large storage requirement combined with present 
computer capabilities severely limits the number of node points that can be used 
and is therefore cause for concern in developing a fully-elliptic computer 
program. In the next section, the finite-difference formulation will be 
simplified to a set of partially elliptic equations. The computer program based 
on a partially elliptic formulation requires much less computer storage than a 
fully elliptic program but only a small increase in storage over that required 
for a parabolic program. 
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III. FOEMULATIONS FOR A PARTIALLY-ELLIPTIC FLOW FIELD 


As discussed in the Introduction, the main flow field in the supersonic 
combustor possesses a predominant flow direction. The diffusion of mass, 
momentum, energy, etc. along the flow direction can be neglected as compared 
with the corresponding convection. The governing equations, except for 
pressure, can be simplified to a set of parabolic equations. By assuming a 
(three-dimensional) pressure field, this set of parabolic flow equations can be 
solved by a marching-integration procedure along the flow direction (assumed to 
be in the z-direction) - The pressure field is generally governed by an 
elliptic equation and can be solved by an iterative procedure. In this section 
the finite-difference equations for such a partially elliptic flow field are 
obtained. The nature of the equations as related to the method of solution is 
discussed. 

Since the flow variables (except pressure) at a streamwise station depend 
only on quantities at the immediately upstream station (not on those at the 
downstream stations), slightly different control volumes are chosen for the 
finite-difference foimulation. The projections of a general control volume 
in the xy and yz planes are shown in figure 3. In the xy ( cros.s-stream) plane, 
the control volume is exactly the same as that in figure 2, whereas, the yz 
(streamwise) plane, the upstream and downstream faces of the control volume 
are located at U and P. By neglecting all diffusion terms in the z-direction in 
equation (l), and taking volume integrations of equation (l) over the respec- 
tive control volumes, a set of finite-difference equations is obtained. The 
continuity equation becomes 


lU 



[(pw)p - (pw)y](Ax)(Ay) + [(pu)^ - (pu)^](Ay)(Az) 


+ [(pv)„ - (pv) ](Ax)(Az) = 0 
n s 


( 9 ) 


The general difrerential equation for 0, except three momentum equations, takes 
the following form. 
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p e w s n 
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( 11 ) 


X X y y 

Where A , A , A‘^ , A , S are the expressions defined in Section II; F depends 

0wsn*u.5(|) U 

on the mass flux through the upstream face of the control volume. 


Py = (pw)^ (Ax)(Ay)/(Az) 


( 12 ) 


For the three momentum equations, the difference equations are fomulated based 
on the control volumes of velocity components 
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where = (Ax) ( Ay ) (Az) /6x^ , = ( Ax ) (Ay ) (Az ) /6y^ , and = (Ax)(Ay). 

When the local Mach nuiriber is greater than one, the last term (pressure gradient) 
in eq. ( 15 ) is replaced hy D^(pp - p^). 

The houndaiy conditions for all flow variables must be specified at the 
boundaries parallel to the z-axis. The formulations of these boundary condi- 
tions were discussed in Section II. At the upstream boundary all flow variables 
must be specified; at the downstream boundai:y, only the pressure is required in 
the subsonic flow field. 

By comparing the finite- difference equations of the parti ally-elliptic 
formulation presented above with these of a parabolic flow problem outlined 
in reference 1 or 13, it is seen that the two formulations are essentially the 
same. The differences lie mainly in the treatment of the pressure field and 
the pressure gradient terms in the three momentum equations. For the parabolic 
flow problem as described in ref. 13, the streamv-ise pressure gradient is 
decoupled from the transverse pressure gradients. Since the downstream pressure 
effects are negligible in parabolic flows, an approximate velocity field is 
obtained by setting dp/dz = 0 or Pp(x,y) = Py(x,y). The approximate pressure 
and velocity fields are then corrected using a pressure correction governed by 


16 


a two-dimensional (elliptic) equation derived from the continuity equation. 

Thus the pressure can he treated as though it were like the other variables, 
independent of downstream quantities. 

In the parti ally-elliptic flow problem the effects of downstream quantities 
transmitted through the pressure field are important. To determine such flow- 
fields , an iterative marching-integration procedure is used. An approxims.te 
three-dimensional pressure field denoted by p* is first assumed and later 
obtained from the previous iteration. The approximate pressure is corrected 
by a pressure correction p' governed by a three-dimensional elliptic equation 
derived from the continuity equation (9)* The finite-difference form of the 
pressure correction equation is 


■^P^P ^s^S ''' ^u 


(16) 


where 


I I t I I t t 

A = A + A + A + A + A_ + A__ 
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= (Ax) (Ay) (pD^/A^)p 


(IT) 
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(I8b) 


(l8c) 


IT 



~ (Ax)(Ay)[(.pD'^/A^)y + (w "1^)^] 


(I8d) 


S = [(pw) - (pw) J(,Ax)(Ay) + [(pu) - (pu) ](Ay)(Az) ( 19 ) 

U U -L' we 

+ [(pv) - (pv) ](Ax)(Az). 

s n 

All variables in eqs. (lT)-(l9) are calculated based on the approximated flow- 
field (p*, etc.). Once the pressure correction p' is obtained from eq. ( 16 ), 
the approximate pressure p* and velocity components (u*, v*, w*) are corrected 
as follows 


p = p* + p' (20) 

u = u* + (d"^/A^)* (p' - p') (21) 

p w p 

V = V* + (D^/A ) (p - p ) (22) 

P s s P 

w = w* + (D^/Ap)* (Pp - p^) + 1^)* Pp (23) 


This procedure is repeated for each forward step. The details of the calcula- 
tion procedures will be described in Section IV. 

By inspection of the above difference equations , it can be seen that the 
solutions at a point (except pressure) depend only on solutions at four 
neighboring points in the xy (cross-stream) plane. The solution in the pressure 
correction depends on quantities at six neighboring points. Therefore, to solve 
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these equations, only the pressure and the pressure correction require three- 
dimensional storage. The other variables require only two-dimensional storage 
like a parabolic flow program. In comparison with the three-dimensional fully- 
elliptic formulation presented in Section II , a great saving of computer storage 
is realized for a three-dimensional computer program based on the partially- 
elliptic fomulation. 
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IV. COMPUTATION PROCEDURES AND COMPUTER PROGRAM 

In this section the numerical procedures for solving the system of three- 
dimensional* parti ally-ellipti e , finite- difference equations are described. 

The computer program thus developed is discussed. 

13 

The present numerical procedure combines the marching-integration scheme 

1 5 

used in the parabolic flow program and the iterative scheme used for solving 
the elliptic pressure field. The general calcuI.ation steps are outlined in 
the following : 

1. The pressure field is first assumed at all three-dimensional storage 
locations . 

2. A marching-integration procedure is initiated which marches through the 

flow field. In equations (l3)-(l5)* the pressure gradient terms are evaluated 

from the assumed pressure field, and the coefficients A, source terms S . said 

u,({) 

Fy(j>y for (|) = u,v,w are evaluated from the upstream conditions. Then the two- 
dimensional difference equations (l3)-(l5) are solved by a double-sweep line-by- 
line iterative scheme using a tri diagonal matrix algorithm. 

3. The newly calculated velocity components are checked for mass conser- 
vation at all node points in the xy (cross-stream) plane. The unbalance of 
mass conservation is necessary for caJ.culation of the pressure correction, 
which is obtained from eq. (l6). Since an iterative scheme is used in solving 
the three-dimensional pressure field, eq. (l6) is simplified by setting 

py = Pp = 0 to save computer storage. Thus, eq. (l6) is reduced to a two- 
dimensional equation, and can be soi.ved by the double-sweep line-by-line itera- 
tive scheme. 
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U. After obtaining the pressure correction, the pressure and velocity 
components are then corrected by eq.s. (20)-(23). 

5 . Equation (lO) for the other variables (e.g. = k, e, H, f, etc.) is 

also solved by the double-sweep line-by-line iterative scheme, so as to pro- 
vide flowfield quantities appropriate to a new downstream station. 

6. The strearawise march is continued until the end of. the partially 
elliptic flow region is reached. By the end of one complete march (iteration), a 
new three-dimensional pressure field has been obtained, 

7 . Steps 2, 3, 4, 5 s 6 are then repeated until the; pressure corrections 
at all node points become smaller than a prescribed value. On the last itera- 
tion, the converged flow variables are printed out. 

As discussed in Section IV, except for pressure, the finite-difference 
formulation for the partially elliptic flow program is essentially the same as 
for the parabolic flow program. In the numerical procedures just described, 
steps 2 to 6 within each iteration are similar to those in the parabolic pro- 
gram. Therefore, the parabolic program provides the main part of the present 
partially elliptic program. The present program considers three-dimensional, 
turbulent, reacting (equllibriimi) flow in a supersonic combustor in the form of 
a rectangular parallelepiped. Each lateral boundary of the parallelepiped can 
be either a solid wall, a symmetry boundary, or a surface with a uniform (or 
free) stream. The wall boundaries are allowed to vary smoothly along the flow 
direction. 

The pressure field in the present program requires three-dimensional 
storage. Thus the strearawise dimension of the pressure field determines, and 
also limits, the maximum number of march steps in each iteration. Pressure 
gradients are evaluated from the assumed pressure field; consideration is also 
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given as to whether the local flow field is subsonic or supersonic. In the 
first iterative solution pressure gradients are assumed to be zero; thus the 
first iterative solution is the parabolic flow solution. The hi^er order 
iterative solutions depend on the pressure fields obtained from the- previous 
iterations. To ensure stability of the iterative procedure a pressure (under) 
relaxation factor is introduced. Since the flow problem is highly nonlinear, 
good criteria to determine a suitable value for this factor have not been 
found. At present this factor is determined by numerical experimentation. 
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V. APPLICATIONS OF THE PROGRAM 


In. this section applications of the present program are presented. Be- 
cause of the lack of good quality experimental data, detailed evaluations of 
the program ty comparing vith experiments are not possible at the present time. 
Therefore, results from the present calculations are compared with results from 
the parabolic program. Such comparisons demonstrate the importance of the 
effects predicted by the partially-elliptic program. Furthermore, previous 
evaluation of the parabolic program has already established its accuracies and 
capabilities; the numerical comparisons- between it and the present program infer 
the usefulness of the present program to supersonic combustor research. 

A numerical example was performed for the flow field of a hydrogen 
Jet mixing with a supersonic airstream in a rectangular duct. As shown in 

figure U, a slightly subsonic (M. = 0.99) hydrogen Jet with urlfomi exit condi- 

J 

tions , w. = 1500 m/sec, T = ^00 K, and p. = 0.101 MPa, is located at the 
J JO 

center of the duct. The ducted airstream at the jet exit plane is also assumed 
to be uniform with w = 990 m/sec, T = 1700 K, and p = 0.101 MPa (i.e. M == 1.2). 

The upper and lower walls diverge slightly to allow for shock and Mach waves to 
develop and interact with the subsonic stream at the center of the duct. Com- 
putations were performed for both reacting and nonreacting cases. 

Streamwise pressure variations along the centerline of the duct 0 and at 
the duct cover C for both partially elliptic and parabolic nonreacting calcu- 
lations are shown in figure 5- The convergent partially elliptic calculations 
were obtained after eleven iterations. Note that streamwise marching steps are 
used on the abscissa in fig. 5 instead, of the usual physical distance; the step 
size used in all present calculations is 0.004 of the duct height. The downstream 
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pressure effects considered by the partially elliptic program are clearly 
demonstrated by the comparison. Because the subsonic region is small, the 
downstream effects change the centerline pressure only about 2 percent. At 
the duct comer, the downstream effects through the boundaiy layer are even 
smaller. Nevertheless, the present partially elliptic program does predict 
correctly the qualitative physical phenomena. 

From the same calculations cross-duct profiles are presented in figures 
6a and b of streamwise velocity w, pressure p, and mass fraction of total 
hydrogen f.^^^ at y = 0 and 0.5 cm along the lines OA and DB in fig. k at step 
32 (z = 0.^68 cm). Because of the small downstream effects, the differences 
between the two calculations are small and mainly occur in the subsonic jet 
region (fig. 6a). As expected, there is no effect in the supersonic flow 
region (fig. 6b). 

The streamwise pressure variations along the centerline of the duct 0 
and the duct comer C for reacting calculations are shown in fig. 7* The 
parti ally-elliptic results are obtained based on 11 iterations. Because of 
the combustion of hydrogen, the centerline pressure increases much more abruptly 
than in the nonreacting case. The large increase in pressure propagates up- 
stream in the parti ally-elliptic calculation and causes as much as lh% of pressure 
difference over the parabolic calculation. Again, the difference in pressures 
at the duct corner is relatively small; however, the pressure variations due 
to the occurrence of shock and Mach waves are clearly sho\m. 

Figure 8 shows the convergence be;havior of the streamwise pressures of 
the same reacting calculation. The streamv.'ise march has been extended to 80 
steps. Results show that the rate of convergence depends on the distance from 
the jet exit. The solution converges faster close to the jet exit than 
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downstream; for example, "beyond step 75 » the centerline pressure has still not 
reached the convergent solution even after 26 iterations. 

The reacting results for cross-stream profiles of the streamv^ise velocity 
w and the mass fraction of total hydrogen , f at y = 0 (line OA in Fig. U) 
and at step 32 are presented in fig. 9a. The profiles of the pressure p and 
the mass fraction of water vapor are shown in fig. 9h. As in the results 
shown in fig. 7, the parti ally-elliptic calculations show differences in these 
profiles over the parabolic calculations , especially the large difference in 
the water vapor distribution. 

Reacting results of cross-stream profiles at y = 0.5 cm (line DB in fig. 4) 
and at step 32 are presented in figs. 10a and b. At this location, the down- 
stream effects predicted by the parti ally- elliptic calculation are much smaller 
than those near the jet centerline. 

Cross-stream profiles at step 52 (z = 0.787 cm) are presented in figures 
11 and 12. Because of the large pressure increase and wide spreading of the 
reacting region, the partially elliptic calculation predicts much larger effects 
than those at step 32. Even away from the centerline, the effects predicted by 
the partially elliptic calculation are both quantitatively and qualitatively 
important . 
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VI. DISCUSSION AND CONCLUSION 


To advance the computational capability for combustor research, a three- 
dimensional partially elliptic flow computer program has been developed. The 
program is capable of predicting combustion flow-fields with large downstream 
effects; this was not possible using a parabolic flow computer program. The 
physical models and numerical schemes of the development of this partially- 
elliptic flow program have been consistent with those of the parabolic flow 
program. 

Physically, this partially elliptic flow program adopts a "k-E" two- 
equation turbulence model and a model of equilibrium reactions among species. 
Although the application of such an equilibrium reaction model is limited, it 
is not difficult to incorporate any finite-rate reaction model into the pro- 
gram. Numerically, this program is based on a finite-difference formulation 
in a "staggered” grid system. A combined central and upwind difference scheme 
is used in the finite-difference equations which are then solved by a Idne-by- 
line iterative scheme using a tridiagcnal matrix algorithm. 

The main difference between the present partially elliptic flow program 
and the parabolic flow program is implied by their names. The parabolic flow 
program is based on a set of parabolic equations with a known (arbitrary) 
pressure field. The partially elliptic flow program is based on a similar set 
of parabolic equations; however, the pressure gradients in the momentum equations 
are evaluated from a three-dimensional pressure solution. Thus an iterative 
marching-integration procedure is necessary to obtain the three-dimensional 
pressure field and the other flow variables . With only a small increase in 
computer storage over. the parabolic flow program, the present program is 
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capable of predicting the downstream effects. 

The importance of the present partially elliptic flow program has been 
demonstrated by the numerical example presented in Section V. The differences 
between results obtained from the partially elliptic and parabolic calculations 
are significant especially in the combustion flow field. Moreover, since 
chemical reactions are strongly pressure dependent, the downstream pressure 
effects predicted may have important effects on calculations using finite-rate 
reactions . 

The main physical limitation of the partially elliptic program is its 
inability to predict recirculation flow fields. For a three-dimensional re- 
circulation flow field, a three-dimensional elliptic flow program is required. 
The finite-difference formulation and method of solution for such an elliptic 
flow have been presented in Section II. However, the large computer storage 
requirement has been a major drawback to development of such a computer program. 

The future extensions for the present program include the optimal applica- 
tion of the iterative procedure in order to minimize computer time. Since the 
nimiber of iterations required depends on the assumed pressure field and the 
relaxation factor, a rational method to estimate a relaxation factor and to 
obtain stable iterations is more important. Moreover, it also has been noted 
that the rate of solution convergence decreases along the streamwise direction: 
i.e. the upstream flow field requires fewer iterations than the downstream. 

Thus , the present program should be tailored especially to actual needs . Other 
work related to combustor development should include the incorporation of a 
■ finite-rate chenn'-cal reaction model to predict realistic combustion flow fields 
and the development of efficient solution algorithms to save computer storage 
and computing time. 
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In conclusion, the present partially elliptic flow program has advanced 
the computational capability for combustor research. Combustor flow fields 
with important downstream effects can now bei predicted. 
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Table I. The Appropriate Exchange Coefficients and Source Terms 

S . for Variable ({) 
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Table II. 

Wall Botindary Conditions 

at Near Wall Points 


^near wall 

^(j) , near wall 

Velocity components 


0 

normal to the wall 

( 0 at wall ) 


Velocity components 
parallel to the wall 

( 0 at wall ) 

+ 

1 + y > 11-5 

~ ln(Ey ) 



y‘^^11.5 

k 

— 

0 


( 0 at wall ) 


e 


- 

f 

- 

0 

H 

- 

n . + 

^ > 11.5 

^t ,H - ZniEy ) + 



p 7" y 

a 


33 




Figure 1.- Grid and control volume for scalar quantities. 
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Partially -elliptic calculation 


Parabolic calculation 


(a) y ^ 0. 

Figure 6.- Cross-duct profiles of velocity, pressure, and mass fraction at step 32 (z=0.46S cm ); 
non-reacting case. 
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(a) Streamwise velocity and total hydrogen mass fraction. 
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(b) pressure and water vapor mass fraction 


Figure 9.- Concluded. 
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Figure 10.- Cross-duct profiles at the 32nd marching step; y = 0.5 cm , reacting case 
(a) Streamwise velocity and total hydrogen mass fraction. 



(b) Pressure and water vapor mass fraction. 
Figure 10.- Concluded. 




Figure 11.- Cross-duct profiles at the 52nd marching step; y = 0. .reacting case, 
(a) Streamwise velocity and total hydrogen mass fraction. 
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(b) Pr6ssur6 and water vapor mass fraction. 
Figure 11.- Concluded. 
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Figure 12.- Cross-duct profiles at the 52nd marching step; y = 0.5 cm , reacting case 


(a) Streamwise velocity and total hydrogen mass fraction. 
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(b) Pressure and water vapor mass fraction. 
Figure 12.- Concluded. 
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